o 
o 

(N 
o 

Q 

(N 



The Polya Urn: Limit Theorems, Polya Divergence, 
Maximum Entropy and Maximum Probabihty 



Marian Grendar 
Department of Mathematics, FPV UMB 
974 01 Banska Bystrica, Slovakia 
Inst, of Mathematics and CS, Banska Bystrica, Slovakia 
Inst, of Measurement Sci., Bratislava, Slovakia 
marian.grendar@savba.sk 



Robert K. Niven 
School of Aerospace, Civil and Mechanical Engineering 
the University of New South Wales at ADFA 
Northcott Drive 
Canben-a, ACT, 2600, Austi-alia 
r.niven @ adfa.edu. au 



o 



> 

O 



X3 
O 

o 



X 



Abstract — Sanov's Theorem and the Conditional Limit Theo- 
rem (CoLT) are established for a multicolor Polya Eggenberger 
urn sampling scheme, giving the Polya divergence and the Polya 
extension to the Maximum Relative Entropy (MaxEnt) method. 
Polya MaxEnt includes the standard MaxEnt as a special case. 
The universality of standard MaxEnt - advocated by an axiomatic 
approach to inference for inverse problems - is challenged, 
in favor of a probabilistic approach based on CoLT and the 
Maximum Probability principle. 

1. Introduction 

Consider an urn containing a; > balls of colors i, 
i = 1, 2, . . . , m; let m be finite. A single ball is drawn from 
the urn, recorded and then returned together with c G Z balls 
of the same color. Assuming — nc < mm{ai,a2, ■ ■ ■ ,am), 
the drawing is repeated n times. This sampling is known 
as the multicolor Polya Eggenberger (PE) urn scheme; c.f. 
[13], [27], [17]. Let = ^ be the relative number of 
times a ball of color i is drawn in n drawings. The vector 
J/" = [i/"^ ^ ly^ ^ . . . ^ i/''^^] will be called type [10], or rt-type 
where necessary to stress that it is induced by n drawings. 
Given the PE scheme, the probability 7r(i/";<7, c) that n-type 
z/" will be drawn is (c.f. [27], [17]): 

iT'^- n™ 1 + c) • • • (ai + (n^ - l)c) 



7r(i/";g,c) 



ml N{N + c)---{N + {n-l)c) 



(1) 



N ' 



where N = ^'iLi ™d vector q consists of qi 
i = l,2,...,m. The Polya Eggenberger (PE) distribution 
(HJ contains, as prominent special cases, the multinomial 
distribution for c = (i.e., random sampling; identically 
and independently distributed (iid) outcomes), the multivariate 
hypergeometric distribution for c = — 1 (i.e., sampling without 
replacement) and the multivariate negative hypergeometric 
distribution for c — 1; c.f. [17]. 

Identify the set of possible colors (outcomes, states) with 
support A:" of a random variable X. Following the notation 
of [10], let V{X) be the set of all probability mass functions 
on X. Let Vn{X) be the set of all possible n-types. Finally, 
let n C V{X) be the feasible set of distributions and n„ ^ 
n n Vn {X). The aim of this work is to examine the Sanov 
Theorem for Polya sampling (i.e., the large deviations behavior 
of 71(1^" G 11; g, c)), its associated Conditional Limit Theorem 



(CoLT) and Gibbs Conditioning Principle (GCP), and connec- 
tions to the Maximum Probability (MaxProb) principle [3], 
[30], [14], [22], [23], [24]. The asymptotic investigations are 
conducted under the assumption that N, f3 = ^ and q may 
change with n in such a way that q{n) q E V{X) and 
/3(n) ^ ^ e (0, 1) as n ^ 00. 

II. Polya divergence 

Let (3 e (0, 1), c^O,p,qe V{X), and q + l3cp > 0. The 
Polya divergence I{p \ \q; /3, c) of p with respect to q is given 
by: 

I{p\\q; I3,c) ^ I{p\\q + I3cp) + ^Iiq\\q + pep) + 



(3c 



l+/3c 
/3c 



log(l + /3c), 



where I{a 1 1 fo) — X^Tii ^^S is the Kullback Leibler (KL) 
divergence [21], with standard conventions. By the conti- 
nuity argument, /(p||g;/3, 0) = I{p\\q). When convenient, 
I{p 1 1 g; /3, c) will be replaced by 1 1 q). 

The following key properties of the Polya divergence are 
needed for later analyses. 

1) Non-negativity. Ip{p\\q) > 0, with equality if and only 
if (iff) p^q. 

2) Lower semicontinuity. Ip{p\ \ q) is lower semicontinuous 
in p, q. If qi > 0, for i = 1, 2, . . . , m, then Polya divergence 
is continuous in the pair p, q. 

3) Convexity in p, q. For any A £ [0, 1] and p,p',q,q', it 
holds that A/|(p||g) + (1 - A)/|(/ Hg') > /|(Ap + (1 - 
X)p'\\Xq+il-X)q'). 

4) Partition inequality. If .A = {Ai, A2, . . . , Ak} is a 
partition of X and paU) = J2teA, P^' iAj) - J2teA, 

j = l,2,...,k, then Ip{p\\q) > Ip{pA 1 1 Ia), with equality 
iff p{i\i G Aj) = q{i\i G Aj), i G Aj, for each j. 

5) Pinsker inequality. If c > 0, the total variation distance 
d{p,q) ^ J2iLi \Pi ~ Qi\ is bounded as follows: I^j{p\\q) > 

Proof: The properties can be established along standard 
lines (c.f. [11], [6], [19]). In particular. Properties 1, 3 and 
4 follow from the log-sum inequality. We outline the proof 
of the Pinsker inequality for the Polya divergence. Since 



the partition inequality holds for Polya divergence, from the 
standard argument (c.f. [6]) it is sufficient to consider the 
binary X with p = [p,l — p] and q = [q,l — q], such 
that p > q, and find out under what restriction on 7 the 
difference g{q) = I^{p\\q)—j(P {p, q) remains negative. Note 
that the difference is for p = q, by Property 1. The first 
derivative g'{q) is g'{q) ^ {q - p){j-^^-^^^^^^-^-^^ + 
■^^^ log ^+0'^^(<i+0'^p) _ 8^1^ Since by assumption q < p, 
q + l3cp < p{l + Pc). If c > 0, 2/ = l + Pc-{q + Pcp) < (1 + 
/3c) (1 — q) and y > 1 — q. Hence, in order to assure negativity 
of the derivative, 87 < j^^^- Setting up 7 = 2{i+0c)^ 
establishes the bound. □ 

III. Sanov Theorem for Polya Sampling 

Topological qualifiers are meant in topology induced on 
the TO-dimensional simplex by the usual topology on M™. 
Following [10], for a set H e V{X) and q € ViX), 
infpen ^ipllq) is denoted by /|(n || g). 

Theorem 1 (Polya Sanov Thm.}: Let 11 be an open set. Let 
q{n) — > q, /3(n) — > /? G (0, 1), as n ^ 00. Then, for n — > 00, 

^log7rKGn;g(n),c) = -/^(n||g). 

Proof: The Method of Types [10] approach to Large 
Deviations will be used. 

For c = 0, the Sanov Theorem is already estabUshed, c.f. 
[26], [10], [7]. The rate function is just the KL divergence, 
i.e.,/«(-|k). 

The case of c will be divided into two subcases: c > 
and c < 0. The following inequalities are needed: 

i) n log n — n< log n! < (n+l) log n — n, valid for n > 6, 

ii) {b - l)log6 - (a - l)loga - (b - a) < logr(6) - 
log r(a) < (6 - i) log 6 - (a - i) log a - (& - a), < a < 6; 
due to [18]. 

For c 9^ A ^ ^ (Z-)™ A f ^ Z-, formula (1) can 
equivalently be expressed as [17]: 



7r(i/";g,c) = 



(2) 



where r(-) is the Gamma function. 

Let c > 0. Note then that the other restrictions under 
which (1) and (2) are equivalent are not active, since —nc < 
min(Q;i, a2, ■ ■ ■ , am)- Applying the inequalities i), ii) to (2), 
i log7r(!^"; g(n), c) is, for n > 6, bounded from above by Un 
and from below by L„: 

rr r/ nil / \ 0/ \ \ (m+l)logn 

Un = -liv"" 1 1 q{n);(3{n), c) + ^ + 



1 1 / 

+ ^(log(l+/3(n)c)--K:iog 

\'i=i 



+ - ( log(l + /3(n)c) 



L„ = -/K||(Z(n);/3(n),c) + 

qi{n) +/3(n)c' 



log 



To establish L„ the standard "trick" of binding log 
from above by — X^^ilogm was used, in addition to i) 
and ii). The non-divergence terms will be denoted u„, 
respectively. 

Let the cardinahty of a set A be denoted |^|. |n| < 
\VniX)\ < (n + 1)™; c.f. [10]. Thus, 

-/(n„||g(n);/3(n),c)+/„ < - log7r(!/" e H; g(n), c) < 

n 

< ""^"^^J"^^^ - /(n„ II qin); f3{n), c) + u„. 

Since m is finite, all terms other than /(•) converge to 

zero as n ^ 00. By assumption, q{n) q, (3{n) — > 
f3 G (0,1). Also, by assumption, qi{n) > 0, for all i, thus 
/(• II •; /3(n), c) is continuous. II is assumed to be open. Thus 
/(n„ II q\ f3{n), c) /(n Wq: f3,c) as n ^ 00. 

For c 7^ A (1 - ^) ^ (Z-)" A (1 - ^) ^ Z", the 
formula (1) can equivalently be expressed as: 



7r(!/";g,c) = 



^! r(l-f 



r 1 



n 

i=l 



r 1- 



Nqi 



(3) 

Let c < 0. Note then that the other restrictions un- 
der which (1) and (3) are equivalent are not active, since 
— nc < mm{ai,a2, . . . ,am)- Applying the inequalities i), 
ii) to (3), the probability i log7r(i^"; q{n), c) can be bounded 
by Un from above and by L„ from below, as L„ = 
-A{u^ \\q{ny,0{n),c)+ln, f/„ = -A{u^ || g(n);/3(n),c) + 
Un, where 

1 1 qin): Pin), c) ^ K log < + ^TY^- 



, , l + /3(n)c 1 
pin)c n 



^ log ( - ^ 



9» 



log 



/3(n)c 



1 



(3{n)c n 



1 

+ - 



^ /?(n)c \ /?(n)c n 



^ g,+/3(n)c;.r / g,+/3(n)ct/f 1 
£^ /3(n)c /3(n)c n 



and In, Un stand for terms that converge to as n — !■ 00. 

Using the same argument as for c > 0, ilog7r(i/" € 
II; qin), c) is bounded 

-AiUnWqin); (3in),c) + In < - log7r(i/" € H; g(n), c) < 

< ""^"^^"^^^ - ^(n„ 1 1 qin); /3(n), c) + u„. 

Since m is finite, the terms other than A( ) converge to 
zero, for n ^ 00. Since A(- 1| g(n);/3(n), c) is continuous, 
the argument used above (case of c > 0) implies that 
A(n„|| ^(n); Pin), c) 7(11 jj g'; /?, c), as n ^ 00. □ 

IV. Polya Conditional Limit Theorem 

The Polya information projection pi/3, c) (Polya I- 
projection, or J|-projection, for short) of g on II is defined 



as p{l3,c) = aiginfp^u Ifj{p\ \ q)- The standard /-projection 
[11] is the special (c — 0)-case of the Polya /-projection. 

The Polya Conditional Limit Theorem (CoLT) is an impor- 
tant consequence of the Polya Sanov Theorem. 

Theorem 2 (Polya CoLT): Let q{n) q, (3{n) {3 E 
(0, 1), as n — » 00. Let 11 be a convex, closed set. Let p{(3, c) be 
the /^-projection of q on 11. Let B{q, e) be the e-ball defined 
by the total variation metric, centered at q. Then for any e > 0, 

lim 7r(z/" e B{p{/3, c), e) | z/" e H; q{n), c) = 1. 

Proof: Let B^{p{/3, c), e) ^ P(A')\/3(p(/3, c), e). Apply 
Polya Sanov Theorem to 7r(i/" e B'^{-)\h'" € n;q(n),c) = 
The decay rate ^B^ {■) \\q) - /^(H || g) > 0. 
Since 11 is, by assumption, convex and closed, by convexity of 
Polya information projection (Property 3) there is unique /^- 
projection of q on 11. Types thus asymptotically concentrate 
on it. □ 
Polya CoLT has the same interpretation as the standard, iid- 
case, CoLT (see [30], [29], [28], [4], [9]): types induced by PE 
sampling, asymptotically conditionally (on the event i^" G 11) 
concentrate on the Polya information projection c) of q 

on n. 

Setting n = V{X) reduces Polya CoLT into its special case: 
the Law of Large Numbers for PE sampling. 

V. Further Results 

By means of the Polya Sanov Theorem and the bounds used 
for its proof, three additional results can be obtained. 

A. Polya Gibbs Conditioning Principle 

For the iid sampling there is a claim, stronger than the 
Conditional Limit Theorem, known as Gibbs Conditioning 
Principle (GCP); c.f. [7], [10], [12]. Alongside of its proof 
[10], the following Gibbs Conditioning Principle for PE sam- 
pling can be established. 

Theorem 3 (Polya GCP): Let q{n) q, l3{n) /3 £ 
(0, 1), as n ^ 00. Let 11 be a convex, closed set. Let p{(3, c) 
be the /^-projection of q on 11. Then for a fixed t, 

lim tt{Xi =Xi,...,Xt=Xt € 11; g(n),c) = 

n— *oo 

t 

1=1 

Loosely put, asymptotically, conditionally upon the event 
j^" £ n, a fixed-length sequence of drawn colors behaves as 
if it was identically and independently drawn from the Polya 
information projection p(f3, c) of q on 11. 

Its f = 1 special case can be established for c = by means 
of the Pythagoras property of the /-projection and the Pinsker 
inequaUty; see [6]. This approach does not carry on to c 7^ 0, 
as the Pythagoras property does not hold for /^-projection 
with c ^ 0. 



B. Maximum Probability - Maximum Entropy Correspondence 

Let v^{fi^c) = argsup,^„gjj^ 7r(!^"; q, c) be the Polya /i- 
projection (/U^-projection, for short) of q on n„; i.e., the 
supremum-probable n-type in n„. Using the J7„, L„ bounds 
(c.f. proof of Polya Sanov Theorem), the asymptotic identity 
of Polya yu-projections and Polya /-projections, can be estab- 
lished along the Unes of [15]. 

Theorem 4 (MaxProb/MaxEnt): Let q{n) q, (3{n) 
P G (0, 1), as n ^ 00. Let Mn{P{n), c) be a set of all I^'f3(^ny 
projections of q{n) on n„. Let X| be a set of all /^-projections 
of q on n. Then, for n ^ 00, 7W„(/5(n), c) = I(/9, c). 

This permits a deeper interpretation of the Polya CoLT. 
Informally: types, asymptotically conditionally (upon z^" G 11) 
concentrate on the most probable type. Even more loosely 
put, the most probable is asymptotically conditionally the only 
possible. 

The asymptotic identity of Polya /i-projections and Polya 
/-projections is illustrated by the following Example. 

Example 1: Let X = {1, 2, 3, 4}, i.e., there are four colors, 
associated with the numbers. Let IV — {p : Yl'i=iPi^i ~ 
3-2,ElliPi = !}• Let n = 10,50,100,1000 and iV(n) = 
100, 500, 1000, 10000, so that /3„ = 0.1, for all considered n. 
Let q{n) = q = [21, 25, 31, 23]/100 for all considered n. 
For each n, let c G {—2, 1, 0, 1, 5, 10}. The Table in Appendix 
A contains the Polya /x-projection of q{n) on n„. In the 
last block of the Table, the Polya /-projection p of 5 on 11 is 
presented, for each considered c. 

C. Polya Conditional Equi-concentration of Types 

The Conditional Equi-concentration of Types (CET) on /- 
projections, an extension of CoLT to the case of 11 admitting 
more than one /-projection, is discussed at [15]. Sinularly, 
CET holds also for Polya /-projections. 

VI. Applications and implications 

Polya CoLT has similar appUcations (and impUcations) as 
the standard iid CoLT (see [9], [15]), but holds in a broader 
context of PE sampling, which encapsulates the iid one. We 
will briefly discuss two applications. 

A. Polya MaxEnt 

The Boltzmann Polya Inverse Problem (BPIP) contains 

the Boltzmann Jaynes Inverse Problem [15] as its special, 
iid, case. BPIP is constituted by the information-pentad 
{X, Nq, n, c, 11} under which the objective is to select a type 
(one or more) from H. Three examples of BPIP are below. 

Example 2: A network containing i = 1,2, ...,m, critical 
nodes, each with ai branches, is accessed by users. Each time 
a connection to node i is accessed, c new connections to this 
node are established (c = — 1 indicates simple congestion, and 
c < — 1, accelerated congestion). Assuming that the number 
of transactions n is known, the objective is to select a type of 
connections to the nodes from 11 = Pn{X). 

Example 3: A stock exchange is estabUshed with A'' 
equiprobable shares, with ai of share i. After a trade in stock 
i, c new shares are issued in it (c < indicates withdrawal of 



shares). The transactions are constrained by the mean value 
of trades in a given period. Given the feasible set n„ of n- 
transactions determined by the mean value of trades, and the 
other above described information, the objective is to select 
an n-type of transactions. 

Example 4: Let n out of N quantum mechanics particles 
be distributed among m energy levels according to the PE 
sampling scheme with initial distribution q and parameter c. 
Let instead of the actual energy distribution (n-type) only 
the mean value of energy of n particles be available. Given 
this information, the objective is to select an n-type from the 
feasible set. 

BPIP is under-determined and in this sense is an ill-posed 
inverse problem. The indeterminacy of the problem translates 
into a multitude of possible methods for its solution. From an 
infinite set of possible methods of solving BPIP, such a method 
has to be selected that, for n — > oo, does not violate Poly a 
CoLT. Clearly, selection of the Polya /-projection of q on 
n„ satisfies the above requirement of asymptotic consistency. 
This selection scheme could reasonably be called the Polya 
Maximum Relative Entropy (MaxEnt) method, where the 
Polya relative entropy is defined as the negative of the Polya 
divergence; H{p 1 1 /3, c) = —I{p \\q;P,c). Note that from 
the point of view of maximization over p, the Polya relative 
entropy effectively reduces to: 

-^Pt logPi +^(pt + ^ logfe + f3cpi). 

The other way for solving/regularizing BPIP that is asymp- 
totically consistent is Maximum Probability (MaxProb), that 
selects the /i^ -projection of q on n„. By Polya Max- 
Prob/MaxEnt the two methods asymptotically coincide, but 
for finite n they make, in general, a different choice. 

The feasible set 11 can for instance (as in the above Exam- 
ples 3 and 4) be formed by moment-consistency constraints 
n = {p : TnLiPiUjixi) = aj,j = 0, 1, 2, . . . , J}, where 
Uj{-) is a given real-valued function, uo(') — 1; o.j is given 
number, ao = 1; such a feasible set is also known as the 
Unear family of distributions. The Polya /-projection of q on 
the linear family of distributions 11 is then impUcitly given by: 

p. (/3, c) = , , , . (4) 

1 ) Distribution of anyons: Not surprisingly, for c = 0, the 
probabiUty distribution (4) turns into the familiar exponential 
(Maxwell-Boltzmann) form of the /-projection on the linear 
family; [II]. For c = —1, the distribution gives the Fermi 
Dirac distribution whilst for c = -|-1, the distribution gives 
the Bose Einstein distribution. These are generalizations of 
the standard Bose Einstein and Fermi Dirac distributions, in 
the sense that a general (not necessarily uniform) sampling 
distribution q is assumed [25]. In this respect it is worth 
recalhng the Example 4 and noting that the PE distribution 
(4) contains an ansatz distribution of quantum-mechanical 
anyons (i.e., particles with properties intermediate between 
those of bosons and fermions; [32]) proposed at [I] as its 



special (uniform q) case. This in our view, provides both a 
probabilistic underpinning of the Acharya & Narayana Swamy 
[1] distribution of anyons as well as its extension to the 
non-uniform sampling case. Further discussion will be given 
elsewhere [25]. 

2) Limitation of axiomatic approach to linear inverse prob- 
lems: We would like to stress that for c ^ the standard 
MaxEnt method [16], [20] (i.e., selection of /-projection 
of q on n„), when applied to BPIP, does not satisfy the 
requirement of asymptotic consistency. Thus, although the 
standard MaxEnt is advocated by an axiomatic approach as the 
logically consistent way of solving ill-posed inverse problems 
with n defined by the moment-consistency constraints (c.f. 
[8], [9]), the method, when applied under PE sampling with 
c ^ 0, violates the Polya Conditional Limit Theorem. This 
reveals a limitation of the axiomatic approach to inference in 
the inverse problems context. 

B. Rare events simulation 

Polya CoLT and Polya GCP can be used for rare events 
simulation in the context of PE sampling, in the same way 
that the standard CoLT and GCP are used in the iid sampUng; 
c.f. [5]. 

VII. Summary 

The standard Conditional Limit Theorem (CoLT) [6] for 
iid sampling provides a probabilistic justification (c.f. [30], 
[9]) of MaxEnt method in the context of so-called Boltzmann 
Jaynes Inverse Problem (BJIP), [15]. In [14] it was suggested 
that MaxEnt can be viewed as an asymptotic instance of 
the MaxProb method, which under the limited information 
available to the BJIP, selects the type (i.e., the empirical 
distribution) with the highest probability of occurrence, from 
the given data-sampling distribution. It was proposed in [30], 
[23], [24], that MaxProb can be considered in a broader 
context; in particular under sampling schemes other than 
the random (i.e., iid) sampling. There it was also pointed 
out that every sampling scheme might be associated with 
its own instance of MaxProb and its own relative entropy 
maximization method. For a particular sampling scheme (or, 
probabiUstic question of certain form, in general) and adjoint 
inverse problem, the relevant entropy maximization can be 
discovered by considering the associated CoLT. The relevant 
CoLT, in turn, provides probabiUstic justification of the asso- 
ciated relative entropy maximization method in the context of 
the inverse problem. Motivated by these observations, in this 
work we have established CoLT for the PE sampling scheme 
and discussed some of its consequences and applications. 

VIII. Notes on literature 

An early physics-motivated work that extends Boltzmann's 
Maximum Probability principle, steps into the direction of 
Sanov Theorem for non-iid sampling and contains a few views 
ahead of its time is Vincze's [30]; see also [31]. For sampling 
without replacement (i.e., (c = — l)-case of PE sampling), the 
Sanov Theorem was established by [12]. A communications 
channel with Polya noise has been considered at [2]. 



TABLE I 

MaxProb to Polya MaxEnt convergence 







n=10 




c=-2 





0.2 


0.4 


0.4 


c=-l 





0.2 


0.4 


0.4 


c=0 





0.2 


0.4 


0.4 


c=l 


0.1 


0.1 


0.3 


0.5 


c=5 


0.1 


0.1 


0.3 


0.5 


c=10 





0.2 


0.4 


0.4 


n=50 




c=-2 


0.06 


0.14 


0.34 


0.46 


c=-l 


0.06 


0.14 


0.34 


0.46 


c=0 


0.06 


0.14 


0.34 


0.46 


c=l 


0.06 


0.14 


0.34 


0.46 


c=5 


0.06 


0.14 


0.34 


0.46 


c=10 


0.06 


0.14 


0.34 


0.46 


n=100 




c=-2 


0.05 


0.15 


0.35 


0.45 


c=-l 


0.06 


0.14 


0.34 


0.46 


c=0 


0.06 


0.14 


0.34 


0.46 


c=l 


0.06 


0.14 


0.34 


0.46 


c=5 


0.07 


0.14 


0.31 


0.48 


c=10 


0.07 


0.14 


0.31 


0.48 


n=1000 




c=-2 


0.056 


0.142 


0.348 


0.454 


c=-l 


0.060 


0.141 


0.338 


0.461 


c=0 


0.062 


0.141 


0.332 


0.465 


c=l 


0.064 


0.141 


0.326 


0.469 


c=5 


0.070 


0.140 


0.310 


0.480 


c=10 


0.073 


0.140 


0.301 


0.486 




P 


c=-2 


0.05628 


0.14179 


0.34759 


0.45434 


c=-l 


0.05974 


0.14126 


0.33826 


0.46074 


c=0 


0.06241 


0.14085 


0.33108 


0.46566 


c=l 


0.06453 


0.14052 


0.32537 


0.46958 


c=5 


0.06997 


0.13968 


0.31072 


0.47962 


c=10 


0.07357 


0.13913 


0.30102 


0.48628 
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IX. Appendix A 

The Table 1 illustrates Maximum Probability to Polya 
Maximum Entropy convergence, for PE sampling; c.f. Sect. 
V.B. 
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